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ABSTRACT 

We investigate shock structure and particle acceleration in relativistic magnetized collisionless pair 
shocks by means of 2.5D and 3D particle-in-cell simulations. We explore a range of inclination angles 
between the pre-shock magnetic field and the shock normal. We find that only magnetic inclinations 
corresponding to "subluminal" shocks, where relativistic particles following the magnetic field can 
escape ahead of the shock, lead to particle acceleration. The downstream spectrum in such shocks 
consists of a relativistic Maxwellian and a high-energy power-law tail with exponential cutoff. For 
increasing magnetic inclination in the subluminal range, the high-energy tail accounts for an increasing 
fraction of particles (from ~ 1% to ~ 2%) and energy (from ~ 4% to ~ 12%). The spectral index of 
the power law increases with angle from —2.8 ± 0.1 to —2.3 ± 0.1. For nearly parallel shocks, particle 
energization mostly proceeds via the Diffusive Shock Acceleration process; the upstream scattering is 
provided by oblique waves which are generated by the high-energy particles that escape upstream. For 
larger subluminal inclinations, Shock-Drift Acceleration is the main acceleration mechanism, and the 
upstream oblique waves regulate injection into the acceleration process. For "superluminal" shocks, 
self-generated shock turbulence is not strong enough to overcome the kinematic constraints, and the 
downstream particle spectrum does not show any significant suprathermal tail. As seen from the 
upstream frame, efficient acceleration in relativistic (Lorentz factor 70 > 5) magnetized {a > 0.03) 
flows exists only for a very small range of magnetic inclination angles (< 34°/7o), so relativistic 
astrophysical pair shocks have to be either nearly parallel or weakly magnetized to generate nonthermal 
particles. These findings place constraints on the models of AGN jets, Pulsar Wind Nebulae and 
Gamma Ray Bursts that invoke particle acceleration in relativistic magnetized shocks. 
Subject headings: acceleration of particles — gamma rays: bursts — shock waves 



1. INTRODUCTION 

Collisionless shocks have been implicated in most non- 
thermal phenomena in the Universe. Models of nonther- 
mal emission from Pulsar Wind Nebulae (PWNe), jets 
from Active Galactic Nuclei (AGN) , Gamma-Ray Bursts 
(GRBs) and supernova remnants (SNRs) require a pop- 
ulation of high-energy particles with power-law spec- 
tra, and collisionless shocks are the primary candidates 
to produce scale-free distributions of energetic particles. 
However, despite decades of research, particle accelera- 
tion in astrophysical shocks is still not understood from 
first principles. 

Shock acceleration of charged particles is thought to 
occur by means of two main mechanisms: Diffusive Shock 
Acceleration (DSA, or first-order Fermi acceleration) and 
Shock-Drift Acceleration (SDA). 1 In first-order Fermi ac- 
celeration, particles stochastically diffuse back and forth 
across the shock front and gain energy by scattering off 
mag netic turbulence embedded i n the converging flows 
(e .g.,|Blandford fc Ostriker|1978| |Bell|1978| |Drury|[l983 



Blandford & Eichler 1987[). SDA is the process whereby 



particles gain energy from the background motional eleo 
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1 In electron-ion magnetized shocks, Shock-Surfing Acceleration 
is also important, a mechanism by which ions repelled by the shock 
electric potenti al are accelerated by the background motional elec- 
tric field (e.g., [Lee et al.|[l9 96). In this work, we consider only 
electron-positron shocks, where the shock-surfing mechanism does 
not operate. 



trie field E = — (3 x B when gyrating across the shock 
front, while they drift parallel (or anti-parallel , accord- 
ing to their charge) to the electric field (e.g., |Chen fc 
Armstrong! |T975] |Webb et aL][T983| |Begelman Kirk 



1990|. SDA is a very fast process compared to DSA 
( |Jokipiij[l982| ); however, it is not very efficient without a 
means for the particles to remain close to the shock and 
encounter it more than once, thereby being repeatedly 
accelerated. It is generally accepted that when a suffi- 
cient amount of magnetic field turbulence exists on both 
sides of the shock, both mechanisms operate in tandem 
making the overall acceleration process very efficient. 

The efficiency of particle acceleration, namely the 
fraction of post-shock particles and energy stored in 
a suprathermal tail of the particle spectrum, may de- 
pend on the composition of the upstream flow (electron- 
positron, electron-ion or pairs-ion plasma), as well as 
its bulk Lorentz factor and magnetization. If the up- 
stream medium is magnetized, an additional parame- 
ter is the obliquity angle that the upstream magnetic 
field makes with the shock direction of propagation. For 
strongly magnetized shocks, cross-field scattering is sup- 
pressed and the particle gyro-centers are constrained to 
slide along the field lines, which are advected downstream 
from the shock. It is then natural to define a critical mag- 
netic obliquity above which particles sliding along the 
magnetic field should be moving fast er than the speed 



of lig ht in order to return upstream (Begelman & Kirk 
19901. Since efficient acceleration requires that particles 



repeatedly cross the shock, such "superluminal" geome- 
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tries should be very inefficient for particle acce leration 
( [Kirk fc Heavens||1989| |Ballard fc Heavens|p9l| , unless 
there is strong scattering across field lines or apprecia- 
ble fluctuations of the magnetic field at the shock, that 
may create local "subluminal" configurations. If the up- 
stream medium is not turbulent by itself, magnetic field 
fluctuations may be generated by the high-energy parti- 
cles accelerated at the shock that escape upstream. In 
turn, such magnetic waves could scatter the particles and 
mediate their acceleration. The highly nonlinear prob- 
lem of wave generation and particle acceleration should 
then be addressed self-consistently. 

Most progress in studying particle acceleration in mag- 
netized shocks has been m ade using semi-analyt i c ki- 
netic theory methods (e.g., Kirk fc Heavens||1989] |Bal- 



lard fc Heavens||T99"Tl |Kirk et al ||2000| |Acht erberg et af 



2001 



Keshet & Waxman 



|2005I 



or Monte Carlo test- 



particle simulations (e.g., |Bednarz fc O strowski 1998 
Niemiec fc Ostrowski|2004]|Ellison fc Double|2004j Botn 
approaches study particle acceleration with simplifying 
(and somewhat arbitrary) assumptions about the nature 
of magnetic turbulence near the shock and the details 
of wave-particle interactions. In particular, it is found 
that superluminal shocks might accelerate as efficiently 
as subluminal configurations, but the required level of 
magnetic turbulence mus t be unreasonably large (Os- 
trowski & Bednarz 2002 1 . The open question is then 
if superluminal shocks can self-consistently produce such 
a high degree of turbulence. 

Fully kinetic particle-in-cell (PIC) simulations provide 
a powerful tool for exploration of the fluid and kinetic 
structure of collisionless shocks on all scales, thus de- 
termining the nature of the shock-generated electromag- 
netic turbulence simultaneously with the spectrum of ac- 
celerated particles. PIC simulations of rela tivistic un- 
maqnetized shocks have been prese nted by Spitkovsky 
(|20'05l |2008a|), [Chang et al.| ( |2008[ ) and |Keshet et at 
( 2008p ; |Spitkov sky ( 2008b]) has shown that unmagnetized 
shocks naturally produce accelerated particles as part of 
the shock evolution. One-dimensional PIC simulations of 
magnetized perpendicular pair shocks, i.e., with the mag- 
netic field along the s hock surface, have shown negligible 
partic le acceleration ( Langdon et al.|[l988 Gallant et al.| 
19921, unless a sizable fraction of the upstrea m bulk kr 



netic energy is carried by heavy ions (Hoshino et al. 1992 
Amato & Arons||2006|) 



In this work, we investigate via first-principles two- 
and three-dimensional PIC simulations the acceleration 
properties of relativistic magnetized pair shocks. We fix 
the upstream magnetic field strength and explore how 
the efficiency and mechanism for particle acceleration 
vary with magnetic obliquity. We find that sublumi- 
nal configurations produce a population of energetic par- 
ticles, and the downstream particle spectrum develops 
a high-energy power-law tail. We confirm that parti- 
cle acceleration in superluminal shocks is extremely in- 
efficient, meaning that self-generated turbulence is not 
strong enough to allow for significant cross-field diffusion 
of the particles. For subluminal shocks, we find that the 
acceleration efficiency increases as the magnetic obliquity 
approaches the boundary between subluminal and super- 
luminal configurations. Also, by studying the trajectory 
of high-energy particles extracted from our simulations, 
we show that at low magnetic obliquities, where the mo- 



tional electric field is still small, the acceleration process 
is mostly controlled by DSA, whereas particles are mainly 
accelerated via SDA when the magnetic field inclination 
is close to the super luminality threshold. 

This work is organized as follows: in ^2] we discuss the 
setup of our simulations and the magnetic field geome- 
try; in <|3]the shock structure and internal physics is in- 
vestigated for some representative magnetic obliquities. 
The main results of our work, concerning the efficiency of 
shock acceleration as a function of magnetic obliquity, are 
presented in fj4] In fj5] we comment on the relative im- 
portance of different acceleration mechanisms (namely, 
DSA and SDA) as a function of the magnetic inclination 
angle. We summarize our findings in fj6]and comment on 
the application of our results to astrophysical scenarios. 

2. SIMULATION SETUP 

We perform numerical simulations using the three- 
dimensional (3D) ele ctromagnetic part icle-in-cell (PIC) 
code TRISTAN-MP ( |Spitkovsky||2005| ) , which is a par- 
allel version of the publicly available code TRISTAN 



(Buneman 19931 that was optimized for studying col- 
lisionless shocks. 

The shock is set up by reflecting a stream of pair 
plasma (the "upstream" plasma) off a conducting wall. 
The interaction between the injected and reflected beams 
triggers the formation of a shock, which moves away 
from the wall. This setup is equivalent to the head-on 
collision of two identical relativistic pair streams, which 
would form a forward and reverse shock and a contact 
discontinuity. In our simulations we follow one of these 
shocks, and replace the contact discontinuity with the 
conducting wall. The simulation is performed in the 
"wall" frame, where the "downstream" plasma behind 
the shock has zero bulk velocity along the direction of 
the upstream flow. We stress that the wall frame of our 
simulations does not in general coincide with the down- 
stream fluid frame, since for oblique magnetic configura- 
tions the downstream plasma may have a velocity com- 
ponent transverse to the upstream flow, as we show in 

m 

We perform simulations in both two- and three- 
dimensional computational domains, and we find that 
most of the shock and acceleration physics is captured 
extremely well by two-dimensional simulations (see ^J6] 
and Appendix |A| for a discussion). Therefore, to follow 
the shock evolution for longer times with fixed computa- 
tional resources, we mainly utilize two-dimensional runs. 
All three components of particle velocities and electro- 
magnetic fields are tracked, however. 

We use a rectangular simulation box in the xy plane, 
with periodic boundary conditions in the y direction 
(Fig. l}. The box is 1024 cells wide (along y), corre- 
sponding to approximately lOOc/wp, where the plasma 
skin depth c/ui p — 10 cells. We also performed limited 
experiments with larger boxes, up to 4096 cells wide, 
obtaining similar results. In the expression above, c is 
the speed of light and tu p = •*/ 4ire 2 n u / (70T1) is the rela- 
tivistic plasma frequency of the upstream flow; n u is the 
upstream number density measured in the wall frame, 70 
the Lorentz factor of the injected plasma, e and m the 
positron charge and mass. We fix c — 0.45 cells/timestep, 
so that the temporal resolution of our simulations is 



0.045 u p 1 per timestep. 



At the final time 9000 ui p 1 of our 
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Fig. 1. — Simulation geometry and configuration of the upstream 
magnetic field. 



simulation runs, the simulation domain is up to ~ 50000 
cells (5000 c/ojp) long along x. 

The incoming pair stream is injected along — x 
with bulk Lorentz factor 70 = 15. The injected 
plasma is cold, with energy spectrum in the proper 
frame given by a three-dimensional Maxwcllian f(j) oc 

7\/7 2 — lexp(— 7/A7) with thermal spread A7 = 10~ 4 . 
Each computational cell is initialized with two electrons 
and two positrons. We also performed limited experi- 
ments with a larger number of particles per cell (up to 8 
per species), obtaining essentially the same results. The 
magnetization of the upstream plasma is a — 0.1, where 
a = B\l (4:TTjomn u c 2 ) is the squared ratio of relativis- 
tic Larmor frequency Q h = eB u /(jomc 2 ) to plasma fre- 
quency uj p for the injected flow (_B U is the upstream mag- 
netic field measured in the wall frame). The upstream 
magnetization is chosen such that perpend i cular shocks 
are not mediated by th e relativistic |Weibel ( 1959 1 insta- 
bility. As discussed by Spitkovsky ( |2005p , the growth of 
this instability is suppresse d in perpendicular shocks if 
a > 10~ 3 ; as we show in £3.1 for parallel shocks the 



Weibel instability may still be important for magnetiza- 
tions as large as a — 0.1. The regime we explore may 
be relevant for PWNe and internal shocks in GRBs and 
AGN jets, as we comment in Sj6] 

Keeping the magnetization fixed, we explore a range 
of magnetic inclinations by varying the angle 9 between 
the shock direction of propagation +x and the upstream 
magnetic field B u , which is initialized to point towards 
the first quadrant of the xz plane, as sketched in Fig. [I] 
(i.e., -B x ,u,-Bz,u > and 9 = arctan[_B ZiU /B XiU ]). The 
magnetic obliquity angle 9 is measured in the wall frame. 
We vary 9 from 6 — 0°, which corresponds to a "parallel" 
shock, with magnetic field aligned with the shock normal, 
up to 9 — 90°, i.e., a "perpendicular" shock, with mag- 
netic field along the shock front. For 9^0°, in the up- 
stream medium we also initialize a motiona l electric field 

e u = -Po x B u alon s -f- Hcrc > Po = -V 1 - V7o * is 

the three-velocity of the injected plasma. 

As a byproduct of the shock evolution, particles and 
electromagnetic waves may propagate upstream from the 
shock at the velocity of light. A fixed simulation box with 
a reflecting boundary on the left and an open boundary 
on the right would not be suitable to study the long-term 
evolution of the shock, since particles and fields would 
be removed from the simulation domain when they cross 
the open boundary, thus artificially suppressing any feed- 
back they may have on the shock. Moreover, in the first 
stages of shock evolution, physical quantities in a large 
portion of the simulation box have not been affected yet 



by the presence of the shock; following their time evo- 
lution would imply a waste of computational resources, 
in terms of both time and memory. Instead, we inject 
particles and fields into the simulation box from a "mov- 
ing injector" , which recedes from the wall (at x wa ii = 0) 
along +X at the speed of light, so that its position at time 
t is roughly Xi n j ~ ct. This ensures that the simulation 
domain consists only of those regions which are in causal 
contact with the initial setup, thus saving on computa- 
tional time. Also, our simulation box expands (via dis- 
crete jumps along +x) whenever the injector approaches 
its right boundary, thus allowing to save memory. In 
summary, a moving injector with an expanding domain 
permits to follow the shock evolution as far as the compu- 
tational resources allow, preserving all the particles and 
waves generated by the shock. 

However, since the shock velocity is smaller than the 
light speed, the distance between the injector and the 
shock increases linearly with time. In all PIC codes, a 
numerical heating instability arises when cold relativistic 
plasma propagates for large distances over the numerical 
grid. The instability is a mixture of grid-Cercnkov in- 
stability - due to the fact that electromagnetic waves 
on a discrete grid have subluminal phase velocity for 
large wavenumbers - and spurious plasma oscillations 
(M. Dieckmann, private communication). We reduce this 
numerical instability by using the fourth-ord er solver for 



Maxwell's equations ( Greenwood et ah] 2004 ) , which im- 
proves the numerical phase speed of the shortest wave- 
modes. To further get rid of the instability, we let our 
moving injector jump backwards (i.e., towards the wall) 
every 10000 time steps, resetting the fields to the right 
of the injector; in the following, we refer to this as a 
"jumping injector" . The first jump happens after 40000 
time steps and we tune the jump according to the shock 
speed so that the injected plasma never propagates far- 
ther than 1000c/cjp from the injector before encounter- 
ing the shock. Although this suppresses the growth of the 
numerical instability - which arises if the injected plasma 
moves farther than ~ 1500 c/lj v from the injector - it is 
not a completely satisfactory solution. In fact, each jump 
may remove from the computational domain some of the 
particles and waves propagating upstream, which would 
not happen if the injector were always receding at the 
speed of light. Also, electromagnetic waves moving away 
from the shock may be reflected back. We have tested 
the effects of a jumping injector by letting the first jump 
happen after 30000, 40000 and 50000 timesteps; in j jO] 
we show that a jumping injector does not significantly 
affect our results. 

3. SHOCK STRUCTURE 

In this section we describe the structure and internal 
physics of relativistic magnetized pair shocks. We fix 
the upstream bulk Lorentz factor (70 = 15) and mag- 
netization (o — 0.1) and vary only the obliquity angle 
9, i.e., the angle between the upstream magnetic field 
and the shock normal, as measured in the wall frame. 
According to the criterion anticipated in |l] a shock is 
"superluminal" if particles cannot escape ahead of the 
shock by sliding along the magnetic field, and they are 
inevitably advected into the downstream medium. The 
critical boundary between subluminal and superluminal 
configurations should satisfy in the upstream fluid frame 



4 



9' clit = arccos(/?g h ), where (3' sh is the upstream-frame 
shock speed, which itself depends on the magnetic obliq- 
uity. This corresponds to 9 cr a — arctanpyo t&xi(9' clit )] in 
the wall frame of our simulations. For 70 = 15 and 
a = 0.1, the shock velocity f3' sh computed from MHD 
equations yields 9 cr - lt ss 34° in the simulation frame, with 
weak dependence on both 70 and a. Among the repre- 
sentative sample of magnetic obliquities that we present 
in this section, 9 = 0° (Figs. [2J and |J, 15° (Figs, g 
and |5) and 30° (Figs. [6] and ^ are subluminal shocks, 
the last one being close to the critical obliquity angle 
6 cxit « 34°; instead, 9 = 45° (Figs. [8] and [9) is a repre- 
sentative superluminal shock. We anticipate (see q6]for 
further discussion) that the main results presented below 
for 70 = 15 and a = 0.1 hold across a larger range of up- 
stream Lorentz factors (we tried 70 = 5 and 70 = 50, for 
fixed a = 0.1) and magnetizations (we tested a = 0.03 
and a = 0.3, for fixed 70 = 15). 

Figs. [2] [4] [6] and [8] present the shock internal structure 
as a function of the longitudinal coordinate x at time 
oj p t = 2250. We include both two-dimensional plots in 
the xy simulation plane and transversely-averaged (along 
y) plots for the following quantities: particle number den- 
sity (panel a for the two-dimensional plot and panel b 
for the corresponding y-average), magnetic energy frac- 
tion e B (defined as the ratio of magnetic energy density 
to kinetic energy density of the incoming flow, panels 
c and d), transverse magnetic field B z (panels e and 
/), transverse electric field E y (panels g and h), longi- 
tudinal electric field E x (panels i and j). We also plot 
the y-avcragcd electric pseudo-potential ~J(E x )dx in 
panel k, normalized to the kinetic energy of injected par- 
ticles. The positron phase space density projected onto 
the x — 7/3*, x — 7/3 y and x — j/3 z planes is shown in pan- 
els I, m, n respectively. Figs. [3][5 [7J[9] show the positron 
7/3 x — 7/3 z momentum space and electron and positron 
energy spectra at three different locations through the 
flow, as marked by arrows at the bottom of Figs . [2j[4j[6j[8] 
respectively. 

In all the simulation runs, a collisionless shock is 
formed which propagates away from the wall. The jump 
in density and electromagnetic fields at the shock, as 
well as the shock velocity, a re in agreement with MHD 
calculations, as we show in |3.1| For oblique configura- 
tions (0° < 9 < 90°, i.e., with the exception of parallel 
and perpendicular shocks), the magnetic field refracts 
at the shock. In the frame of the simulation, the lon- 
gitudinal speed of the incoming flow drops to zero af- 
ter the shock; however, for oblique configurations, the 
downstream fluid preserves a small transverse velocity 
Pz.d < 0, in agreement with MHD predictions. In this re- 
spect, the wall frame of our simulations does not strictly 
coincide with the downstream fluid frame, but for most 
purposes /3 Z)C j is small enough that the two frames are 
almost undistinguishable. In the shock transition layer, 
the incoming cold beam isotropizes and thermalizes to 
a downstream temperature that relates to the upstream 
kinetic energy, reduced by the fraction of upstream par- 
ticle energy converted to downstream magnetic energy. 
For low obliquities (9 < 45°), the downstream plasma 
is fully isotropic in three dimensions. For larger obliq- 
uities, the plane perpendicular to the downstream mag- 
netic field, where particle orbits would lie in the absence 
of significant motion along the field, is almost degenerate 



with the direction of propagation of the incoming plasma. 
In this case, particle motions will mostly be confined to 
the plane orthogonal to the downstream field, which pre- 
vents efficient isotropization along the direction of the 
field. In other words, for low obliquities the effective adi- 
abatic index of the downstream plasma is appropriate for 
a three-dimensional relativistic gas (T = 4/3), whereas 
for 9 > 45° the adiabatic index tends to that of a two- 
dimensional relativistic fluid (T — 3/2). 

This section is structured as follows. In §3.1| we present 
a general description of the structure of oblique shocks. 
We discuss that Weibel-like filamentation instabilities are 
likely to mediate the shock formation for low magnetic 
obliquities, whereas magnetic reflection from the shock- 
compressed magnetic field triggers high-obliquity shocks. 
For 0° < 9 < 90°, the longitudinal profile of particle 
density, magnetic energy and transverse magnetic field 
shows a fluid structure with two shocks, a strong "fast" 
shock (which we identify with the main shock) and a 
weak "slow" shock. The velocity and jump conditions of 
the fa st sh ock are i n agreement with MHD calculations. 

In §3.2| and ^ 3.3 we describe in detail some represen- 
tative examples of subluminal (0 < 9 < cr ;t) and super- 
luminal (0 C rit ;$ < 90°) shocks. In subluminal config- 
urations, some shock-accelerated particles escape ahead 
of the shock, and they populate a diffuse beam of high- 
energy particles propagating upstream. These "return- 
ing" particles trigger the generation of upstream waves, 
which propagate towards the shock with a wavevector 
oblique to the background magnetic field. In superlu- 
minal shocks, no returning high-energy particles or up- 
stream waves are observed. The importance of these 
oblique modes for particle acceleration will be discussed 
in fj5j where we present a detailed investigation of the 
acceleration mechanisms for subluminal shocks. 

3.1. Oblique shocks: general overview 

Unmagnctizcd relativistic shocks are mediated by 
Weibel-like filamentation instabilities ( Weibel 11959 



Medvedev fc Loeb|1999||Gruzinov fc Waxman|1999f 7Tn 

coming particles are scattered by the Weibel- generated 
fields and the incom ing flow is isotro pized, decelerated 
and compressed (e.g., Spitkovsky|2005 1 . The same mech- 
anism is likely to mediate moderately magnetized (<r < 1) 
parallel shocks. In fact, we find that the width of the 
transition layer in a = 0.1 parallel shocks, as evaluated 
from the density profile (Fig. |2]& at x B h ~ 730c/w p ) 

)f skin 



of the order of a few tens of skTn depths, as observed in 
unmagnetized shocks. 2 Also, the magnetic energy den- 
sity peak at the shock (reaching ~ 30% of the upstream 
kinetic energy density, see Fig. [2]rf) is a characteristic 
signature of magnetic field amplification by a Weibel-like 
process. The field then decays to the plateau predicted 
by MHD calculations on a length scale ~ 50c/w p . 

Oblique shocks (9 7^ 0°) may be triggered by an- 
other process, namely magnetic reflection of the i ncom- 
ing flow from t he shock-compressed magnetic field ( |Alsop| 
fc Arons||1988l ), as reported fo r perpendicular shock's by 
previous one-di mensional (e.g., |Langdon et al.|1988||Gal-| 



2 We point out that the shock transition region in strongly mag- 
netized (a > 1) plasmas is much wider, since cross-field motion is 
suppressed and the Weibel instability is quenched. 
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PIC simulations. The synchrotron maser instability, ex- 
cited by coherent gyration of the incoming particles in 
the shock magnetic fie ld, mediates the thermal ization of 
the upstream plasma (Hoshino & Arons 1991). Indeed, 
the peak in the transversely-averaged density profile ob- 
served for 8 = 45° at the shock location (see Fig. [8]& 
at x s h w 920c/oj p ) is due to Larmor gyration of the in- 
coming particles in the compressed shock fields, which 
decelerates the upstream flow and therefore increases its 
density. The density rise is as sharp as a few Larmor 
radii in the compressed field (i.e., a few plasma skin 
depths) , much thinner than for parallel shocks (compare 
Fig. [8]& with Fig.[2]&). In just a few Larmor scales down- 
stream from the shock, the density profile flattens to a 
plateau in reasonable agreement with MHD jump condi- 
tions. Shocks with larger obliquities (8 > 45°) arc also 
triggered by magnetic reflection. 

Both Weibel instabilities and magnetic reflection play 
a role in mediating shocks with intermediate obliquities 
(0° < 8 < 45°). Filamentation instabilities dominate 
for lower obliquities, such as 8 = 15°, whose fluid struc- 
ture resembles 8 = 0° in the magnetic overshoot at the 
shock (Fig.|4]rf at ir sll w 720c/cj p ). Instead, such Weibel- 
induced magnetic energy excess is not seen for 8 = 30° 
(Fig. [6]d at x s h » 770c/w p ) or larger obliquities; rather, 
investigation of early times (uj p t = 225) for 8 = 30° shows 
that this shock started out with the prominent density 
peak characteristic of highly-oblique shocks, which sug- 
gests that this shock is mainly triggered by magnetic re- 
flection. For 8 — 30°, the sharp density peak observed at 
early times was later smoothed by the increased pressure 
in shock-accelerated particles. 

For all oblique magnetic configurations (0° < 8 < 90°), 
a fluid structure with two shocks (a strong "fast" shock 
and a weak "slow" shock) is seen in our PIC simula- 
tions. The presence of two shocks is predicted fo r the 
parameters (70 = 15 and a — 0.1) that we adopt (Lan- 



dau et al. 19841 and has been confirmed by MHD sim- 
ulations of the collision of a magnetized shell against a 
wall (S. Komissarov, private communication). The par- 
ticle number density rises across both shocks; however, 
while the fast shock involves an amplification of the mag- 
netic field, across the slow shock the magnetic energy de- 
creases. The slow shock propagates into the downstream 
medium of the fast shock. Its transition layer in our 
PIC simulations is not as narrow as MHD would predict, 
due to insufficient dissipation along the magnetic field. 
Rather, the slow-shock transition region stretches self- 
similarly with time; the jump in fluid quantities between 
the wall and the fast-shock downstream stays roughly 
constant in time. We observe that the slow shock is 
slower and weaker as 8 approaches either 0° or 90°, and 
relatively stronger and faster for intermediate obliquities 
(8 w 30°). For oblique magnetic configurations, we shall 
refer to the fast shock simply as "the shock" and to its 
downstream region as "the downstream medium" . 

For the main (i.e., fast) shock, the velocity and jump 
conditions extracted from our PIC simulations agree 
within a few p ercent with the solution of relativistic 
MHD equations (Appl & Camenzind 1988) supplemented 
by the three-dimensional Syngc ( 195?]) equation of state, 
provided that MHD results are Lorentz-transformed from 



the standard shock frame into the simulation frame. 3 
In Table [JJ we compare the results from PIC simula- 
tions and MHD calculations for the shock velocity /3 s h 
and the jump across the shock in particle number den- 
sity (na/n u ), magnetic energy density (B^/B^), trans- 
verse magnetic field (B z ^/B z u ) and transverse electric 
field (Ey t d/Ey tU ). We remark that, for oblique shocks 
(0° < 8 < 90°), the background motional electric field 
(last column in Table [lj does not drop to zero after the 
shock since, as discussed above, the downstream flow pre- 
serves a small transverse velocity /3 Zj( j < which gives a 
motional field E y & = — [3 Z ^B x d > 0. 

As shown in Table [lj for 8 < 30° our kinetic calcula- 
tions are in remarkable agreement with three-dimensional 
MHD results, suggesting that for low magnetic obliq- 
uities the downstream plasma in our two-dimensional 
simulation box has a three-dimensional effective adia- 
batic index, as anticipated above. For 8 > 45°, the 
adiabatic index of the downstream fluid in our simu- 
lations tends to that of a two-dimensional gas, which 
explains the discrepancy between PIC simulations and 
three-dimensional MHD calculations for 8 = 45° in Ta- 
ble [Q 4 



TABLE 1 

Comparison between PIC simulations (upper row) and 
MHD calculations (lower row) 
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Another common feature of all oblique and perpen- 
dicular shocks (0° < 8 < 90°) is the presence of a 
transient electromagnetic precursor (seen in B z and E y 
for the magnetic configuration shown in Fig. IT]) propa- 
gating at the speed of light ahead of the shock, as re- 



ported in previous one-dimensional (e.g ., |Gallant et al. 
|1992||Hoshino et al.|1992[|Hoshino| 2008| ) and short three- 



dimensional ( Spitkovsky 2005[) PIC simulations of per- 
pendicular shocks. This is a transverse electromagnetic 
wave which is generated in the initial stages of shock 
evolution by bunching of coherent Larmor orbits near 
the shock front. In multi-dimensional simulations, we 
observe that it is no longer produced as soon as the co- 
herence between different locations along the shock sur- 
face is lost. The transient precursor generated at early 

3 Shock speed and jump conditions from our PIC simulations 
are computed at relatively early times (oj p i = 2250), since for later 
stages of the shock evolution the pressure by accelerated particles 
may significantly affect the shock structure. 

4 We point out that the downstream adiabatic index in high- 
obliquity shocks tends to that of a two-dimensional fluid also in 
three-dimensional PIC simulations, since particle isotropization 
along the field is still harder than orthogonal to the field. How- 
ever, in our two-dimensional runs, the simulation plane for high 
obliquities is almost degenerate with the plane orthogonal to the 
downstream field, where particle gyro-orbits lie, and this makes 
such isotropization even harder. 
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times survives as a dispersive wave packet propagating 
upstream, farther from the shock than Figs. [4] [6j [8] show. 
We find that its power increases with magnetic obliquity 
and its Poynting flux may be comparable or in excess of 



3.2. Subluminal shocks: < 6 < 6 C 



t = 2250 of shocks 



The internal structure at time 
with magnetic obliquities 9 = 6° (Figs. [2][3]), 15° 
(Figs.[i][5| and 30° (Figs. ^§f§ shows some features com- 
mon to all subluminal magnetic configurations. Sub- 
luminal shocks are characterized by the presence of a 
diffuse stream of high-energy shock-accelerated particles 
that propagate ahead of the shock and trigger the growth 
of oblique waves in the upstream medium. In turn, as we 
explain in Sj5] such waves are of fundamental importance 
for particle acceleration. 

The longitudinal phase space plots of positrons in pan- 
els l-n of Figs. |2]|4][6]show the injected plasma as a cold 
dense beam propagating with 7/3 x ~ —15 and broad- 
ening while it approaches the shock, which is located 
at x s h pa 730c/o7p for 6 = 0°, at x s h ~ 720c/o7 p for 
9 = 15° and at x sh « 770 c/lu p for 9 = 30°. In the 
shock transition layer, the incoming plasma is isotropized 
and thermalized. A population of high-energy particles 
moving ahead of the shock appears in the upstream re- 
gion as a hot diffuse beam with 7/3 x > (panel / in 
Figs. [2]|4j[6|; as we show in <|5] these particles have been 
acceleratea at the shock. Additional insight is provided 
by the positron 7/3 x — 7/3 z momentum space and elec- 
tron (blue) and positron (red) spectra in Figs. plplR ex- 
tracted from three different locations across thenow as 
marked by arrows at the bottom of Figs. [2] [4] [6] respec- 
tively. Upstream slabs (panels c and / in Figs. [3][5j|7| 
show two components: the cold incoming plasma ap- 
pears in particle spectra (panel / in Figs. p][5j[7| as a 
relatively narrow bump peaking at low energies, around 
7 rj 15, which broadens as the flow approaches the shock; 
the returning diffuse beam of shock-accelerated particles 
populates the high-energy bump, which is more and more 
depleted (especially at intermediate energies) as we move 
farther upstream. Far downstream slabs (panels a and d 
in Figs, mum show roughly isotropic particle distribu- 
tions with thermal spectra supplemented by high-energy 
tails. For 9 — 30°, the downstream population of ener- 
getic particles significantly varies along x 7 with regions 
closer to the shock showing more particles of higher ener- 
gies than slices farther downstream (Fig. 6 l-n, and com- 
pare panels a, d with panels b, e in Fig. 7 1. Since fluid 
elements which are farther downstream passed through 
the shock at earlier times, this suggests that the shock 
acceleration efficiency increased over time, as we discuss 
ing2] 

For oblique magnetic configurations (Fig.Hfor 9 — 15° 
and Fig. 16] for 9 = 30°), the phase space plots x — j(3 y 
(panel m) and x — 7/3 z (panel n) are not symmetric 
around 7$ = (i = y, z), when focusing on the high- 
energy particles located in the vicinity of the shock (see 
also panels b, c in Figs.p][7|. The high-energy positrons 
with 7/3 z > in the snock and upstream regions are 
shock-accelerated particles whose gyro-center is heading 
upstream along the oblique magnetic field (we remind 



that B X ,B Z > 0). 5 With increasing magnetic obliquity 
from 9 = 15° to 30°, the asymmetry around 7/3 z = 
becomes more evident and the returning high-energy 
positrons show increasingly larger t/3 z than 7/3 x (com- 
pare panels I and n between Fig. [4] and Fig. [6]). The 
electron phase space x — 7/3 z (not plotted) shows the 
same asymmetry. 

The excess of high-energy positrons with 7/3 y < ob- 
served in the vicinity of the shock (panel m of Figs. [3] and 
R3I) is a natural consequence of the \7B drift experienced 
by particles accelerated at the shock. For positively- 
charged particles, the drift is directed along — y, which 
explains the sign of the asymmetry in Figs. |4]m,[6]m. 
Since the magnetic field jump at the shock is larger for 
increasing obliquities (Table IT]), and the drift velocity is 
then higher, we expect the asymmetry to be more signif- 
icant for 9 — 30° than for 15°, as observed. Moreover, 
since electrons and positrons move in opposite directions 
under V B drift, the electron x — 7/3 y phase space (not 
shown) exhibits a reversed asymmetry, with an excess of 
particles having 7/3 y > 0. 

The returning beam of high-energy particles triggers 
in the upstream medium the generation of waves whose 
wavevector is oblique to the background magnetic field. 
Such waves have an electromagnetic component (B z in 
panel e and E y in panel g of Figs. [2]|4][6]) and a smaller 
electrostatic component (E x in paneli of Figs. [2j[JJ[6]). 
The electrostatic component also appears as upstream 
oscillations in the transversely-averaged (i? x ) (panel j 
of Figs. 2JEJI6]) and in the electric pseudo-potential 
— J (E x ) ax (panel k of Figs. l2]w]k3]). The upstream waves 
grow in amplitude while propagating towards the shock, 
where the regular oblique pattern observed far upstream 
(with characteristic wavelength ~ 20 c/u p at time uj p t = 
2250) changes to a more complicated structure (e.g., 
compare x > 1000 >c/u} p with 750c/w p < x < 900c/cj p 
for 9 = 15° in Fig. He, g, i). The ratio between the mag- 
netic energy of the oblique waves (mostly contributed 
by their B z ) and of the background upstream field is of 
the order of a few tens of percent, roughly constant with 
magnetic obliquity. 

By comparing panels e, g, i with panels l-n of Fig.|6]for 
9 = 30°, the growth of the upstream oblique modes ap- 
pears to be spatially coincident (at x w 1000 c/u> p ) with 
the leading edge of the population of returning parti- 
cles. The same holds true for smaller obliquities, but for 
9 = 0° and 15° the returning particles that are farthest 
upstream are outside the x-range plotted in Figs. [2] and [4] 
This suggests that interplay between the unperturbed in- 
coming plasma and the returning high-energy beam may 
be responsible for triggering the oblique waves. 

The growth of the upstream waves is neither an un- 
physical consequence of our jumping injector nor an ef- 
fect of the numerical heating associated with the prop- 
agation of a cold beam through the grid. In fact, if we 
replace the reflecting wall at x wa ii — with an open 
boundary, everything else being equal, we observe that 

5 For 8 = 30°, the x- velocity of a particle sliding at the speed 
of light along the upstream magnetic field would be « 0.5 c in the 
simulation frame. This is compatible with the longitudinal location 
x Ri 1100 c/uip at time u) p t = 2250 for the returning particles which 
are farthest upstream (Fig. Ifw-ra). The same is true for smaller 
obliquities, but the leading edge of the population of returning 
particles is outside the x-range plotted in Figs. [2] and [I] 
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Fig. 2. — Shock structure and positron phase space at time u> p t = 2250 for = 0°: a) Particle number density in the simulation plane, 
normalized to the upstream value; 6) Transversely-averaged number density; c) Magnetic energy fraction eb, i.e., the ratio of magnetic 
energy density to the kinetic energy density of the incoming flow; d) Transversely-averaged (eb); e)-/) B z and transversely-averaged (Bz), in 
units of the upstream background magnetic field; g)-h) E y and transversely-averaged (E y ), in units of the upstream background magnetic 
field; E x and transversely-averaged {E x ), in units of the upstream background magnetic field; k) Electric pseudo-potential — f {E x ) dx, 
normalized to the kinetic energy of incoming particles; l)-n) Longitudinal phase space plots x — 7ft (i = x, y, z) of positrons, with phase 
space density plotted as a two-dimensional histogram. For the two-dimensional plots of B z , E y and E x , q a stands for (9/|<?|)|<7| a . 




Fig. 3. — Positron momentum space and electron (blue) and positron (red) spectra at time uipt = 2250 for 9 = 0°, at three slices through 
the flow as marked by arrows at the bottom of Fig. [2] a)-c) Momentum space 7/3 x — -y/3 z in the three slabs, with momentum space density 
plotted as a two-dimensional histogram; d)-f) Electron (blue) and positron (red) spectra 7d7V/d7 in the three slices, normalized to the 
total number of particles within the leftmost slab (apart from a constant multiplication factor). 
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Fig. 4. — Shock structure and positron phase space at time oj p t = 2250 for 9 = 15°. The fluid quantities are normalized as in Fig. [2] 
but here B z and E y are in units of the transverse component £? z , u of the upstream background magnetic field. The excess of magnetic 
energy at the shock (panel d at x s \y 720c/u) p ) does not show up in the transversely-averaged {B z ) (panel /) due to cancellation between 
magnetic fluctuations of opposite sign. 
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Fig. 7. — Positron momentum space and electron (blue) and positron (red) spectra at time ui p t = 2250 for 8 = 30°, at three slices 
through the flow as marked by arrows at the bottom of Fig. [6] See the caption of Fig. [3] for details. 
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no upstream waves are produced, nor is there a shock. 
In order to test the role played by the returning high- 
energy particles in generating the upstream waves, we 
have performed a suite of simulations in which we arti- 
ficially suppress particle acceleration. In these runs, we 
introduced cooling, so that particles with 7 > 7 coo i lose 
a random fraction of t heir nonthermal ene rgy (a similar 
technique was used in Keshet et al. 2008 1 . By varying 
7cooi ( w e tried 7 coo i = 50, 80, 140 for 9 = 0°), we con- 
clude that particle cooling suppresses the growth of the 
upstream waves, and by decreasing 7 coo i (i.e., by forc- 
ing more particles to cool) both the wave amplitude and 
the distance from the shock where they are triggered de- 
crease. 

In summary, the oblique waves seen in the upstream 
medium are generated by the returning beam of high- 
energy particles; they may result from a mixed-mode 
streaming instability which couples the electrostatic two- 
stream instability with the electroma gnetic filamentation 
(Weib el) instability, as envisioned by |Bret fc Dieckmann 
( 2008 1 . A complete analysis of the wave dispersion re- 



lation and generation mechanism is deferred to a future 
paper. 

Finally, we point out that the electrostatic compo- 
nent of the upstream waves induces oscillations in the 
7/? x -velocity of the incoming plasma, as seen in the 
positron phase space of Fig. |2]Z. The electron x — 7/3 x 
phase space (not shown) presents oscillations of oppo- 
site sign. The difference at low energies between electron 
and positron spectra in the upstream medium (panel / in 
Figs. [3]|5]|7| can be correlated with the sign of the elec- 
tric pseudo-potential at the same longitudinal location 
(panel k in Figs. [2j|4]|6] respectively), since a net electric 
potential accelerates incoming charges of one sign and 
decelerate charges of the opposite sign. We further com- 
ment on the effects on particle spectra of t he ele ctrostatic 
component of upstream oblique waves in |4.2| 



3.3. Superluminal shocks: cl -i t < 9 < 90° 

Figs. M and [9] show the shock structure, phase space 
plots and particle spectra for a representative superlumi- 
nal shock (9 = 45°) at time oj p t — 2250. No significant 
particle acceleration occurs in superluminal shocks, and 
the upstream plasma is not appreciably perturbed by the 
presence of the shock. 

In the phase space plot of Fig. [8]/, the incoming plasma 
is seen as a cold beam with 7/3 x « —15, which is 
slightly heated while approaching the shock and then 
sharply thermalized in the shock transition layer at x s h ~ 
920c/w p . As d iscu ssed bylHoshino & Arons (19911 and 
anticipated in §3.1| the incoming plasma thermalizes due 
to the absorption of cyclotron waves produced by co- 
herent particle gyration in the compressed shock fields. 
Thermalization proceeds with time in the downstream 
medium: panels I and m of Fig. [8] show that the plasma 
farther downstream, which had more time to thermalize 
after passing through the shock, has a broader x — j{3 x 
and x—"fP y distribution than the fluid closer to the shock. 

Downstream particles preferentially move in a plane 
perpendicular to the downstream oblique field; indeed, 
the oblique boxy shape of their downstream 7/3 x — 
7/9 z momentum space (Figs. [9]a, b for two different 
downstream slabs) is suggestive of particle gyro-motion 
around the oblique field. As discussed at the beginning 



of ££3j for magnetic obliquities 9 > 45° it is harder for 
downstream particles to isotropize along the magnetic 
field rather than in the plane orthogonal to it. Since 
for high obliquities such plane approaches the xy simula- 
tion plane, this should result in the downstream thermal 
width of the x — 7/3 z phase space being smaller than the 
x — 7/3 x or x — 7/3 y distributions, as observed (compare 
panel n with panels I and m in Fig.M). 

The net fluid velocity /3 Z ^ < behind the shock in 
Fig. [8]n is responsible for the residual motional electric 
field E y ^ > observed in the simulation frame down- 
stream from the shock (hard to see in Fig. but sec 
Table [I]), a feature common to all oblique magnetic con- 
figurations, as discussed in £j3J This residual fluid ve- 
locity vanishes farther downstream, in regions already 
overtaken by the slow shock. The slow shock appears as 
a smooth increase in density (Fig. [8]&) and decrease in 
magnetic energy (Fig. [8]e~) and B z (Fig. [8jf ) behind the 
main shock towards the wall. Similar properties for the 
slow shock are observed for 9 = 30° in Fig. [6j whereas 
the slow shock for 9 — 15° is outside the cc-range shown 
in Fig. |4j 

Most notably, no evidence for particle acceleration is 
present for 9 = 45°. No high-energy particles are ob- 
served to propagate upstream from the shock (Fig. J8]'- 
n), and the particle energy spectrum in upstream slices 
is just the thermal distribution of the injected plasma 
Lorentz-boosted along — x (Fig. [91/* ). The range of down- 
stream 4-velocities is much smaller than for subluminal 
shocks (compare panels l-n in Fig. |8] with the same pan- 
els in Figs. I2]|4]l6]) , and the downstream energy spectrum 
(Fig.^d, e) resembles a purely thermal distribution with- 
out any suprathermal tail. The lack of returning particles 
explains why no oblique waves are seen in the upstream 
medium of superluminal shocks (compare panels e, g, i 
in Fig. [8] with the same panels in Figs. [2]|4][6j . 

We remark that the main features discussed here for 
9 = 45°, and most importantly the absence of shock- 
accelerated particles, are common to all superluminal 
configurations up to the perpendicular case 9 = 90° in- 
vestigated via PIC simulatio ns by, e.g., |Gallant et al.| 
( p92| and |Spitkovs"ky] ( [2005] . 



4. PARTICLE ENERGY SPECTRA 

We now investigate particle acceleration in subluminal 
and superluminal shocks by studying the downstream 
energy spectra of our representative sample of magnetic 
obliquities: 9 — 0°, 15°, 30° (subluminal shocks) and 9 = 
45° (superluminal shock) . We remind that the boundary 
between subluminal and superluminal configurations is 
#crit ~ 34 °. 

In {4.1 we show that subluminal shocks (0 < 9 < # cr it) 
produce a population of high-energy particles as part of 
the shock evolution. Such shock-accelerated particles 
populate a suprathermal power-law tail in the energy 
spectrum behind the shock. As the magnetic obliquity 
increases within the subluminal range, the high-energy 
tail becomes flatter and it accounts for an increasing frac- 
tion of downstream particle number and energy. In other 
words, the acceleration efficiency increases with magnetic 
obliquity for 9 < # C rit- Instead, the downstream parti- 
cle spectrum for superluminal shocks, 6> cr i t < 9 < 90° 
(including perpendicular shocks), does not show any sig- 
nificant suprathermal tail, meaning that particle acceler- 
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Fig. 9. — Positron momentum space and electron (blue) and positron (red) spectra at time uj p t = 2250 for 6 = 45°, at three slices 
through the flow as marked by arrows at the bottom of Fig. [8] See the caption of Fig. [3] for details. 
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ation is extremely inefficient for superluminal configura- 
tions. 

In §4.2| wc follow the time evolution of downstream par- 
ticle spectra and show that the slope and normalization 
of the high-energy tail in high-obliquity - yet still sublu- 
minal - shocks varies with time much more than for low 
obliquities. This suggests that a different acceleration 
mechanism is operating in shocks close to the superlumi- 
nality threshold as compared to quasi-parallel shocks, as 
we discuss in ^J5j 

In the following, all the particle spectra are com- 
puted in a slab 200 c/cj p -wide which is centered 400 c/w p 
downstream from the shock, regardless of time or mag- 
netic obliquity. Strictly speaking, our particle spectra 
are not computed in the downstream fluid frame, due 
to the net fluid velocity P z ^ observed in the simula- 
tion frame behind oblique shocks; however, |/3 z ,d| < 0.1, 
so that our wall-frame particle spectra are an excellent 
proxy for their downstream-frame counterparts. In or- 
der to estimate the acceleration efficiency of sublumi- 
nal shocks, we fit the downstream particle spectrum 
with a three-dimensional Maxwellian (as appropriate for 
8 < 45° shocks, see Q with spread A^ MB plus a power 
law starting at 7 m i n with exponential cutoff at 7 cut 
whose width is A7 cut . The fitting function is /(7) = 
Nus'fy/j 2 - lcxp(-7/A7 MB ) + A^ PL 7 _p min[l, cxp(-(7 - 
7cut)/A7 C ut)], with N PL = for 7 < 7, nin . The high- 
energy tail is defined to start at the intersection of the 
fitting Maxwellian and the fitting power law. 

4.1. Dependence on magnetic obliquity 

Fig. [TO] shows the main result of our work, namely 
the dependence of the acceleration efficiency on mag- 
netic obliquity: downstream spectra for 8 = 0° (blue), 
15° (green), 30° (red) and 45° (black) are shown at the 
final time Lu p t = 9000 of our simulation runs. Subpanel a 
shows the slope p of the fitting power law as a function 
of magnetic obliquity, while subpanels b and c present 
how the fraction of particles and energy stored in the 
suprathermal tail depends on the magnetic inclination. 

For subluminal shocks, the downstream energy spec- 
trum consists of a low-energy thermal bump and a high- 
energy suprathermal tail. The low-energy part, which 
accounts for most of the energy, is populated by the par- 
ticles that are directly transmitted downstream, where 
they are isotropized and thermalizcd. A fraction of the 
incoming particles may be scattered or reflected back up- 
stream when they first encounter the shock; they remain 
in the vicinity of the shock for longer times and may be 
significantly accelerated (see the high-energy population 
at the shock in panels l-n of Figs. |2]|4][6l. If they fi- 
nally escape upstream, they may contribute to the beam 
of returning particles which trigger the generation of 
the upstream oblique waves; if transmitted downstream, 
they populate the high-energy suprathermal tail of down- 
stream particle spectra. 

From Fig. [lO] it is apparent that the suprathermal tail 
in subluminal shocks (blue, green and red lines for 8 = 0°, 
15° and 30° respectively) is much more pronounced than 
in superluminal shocks (black line for 8 = 45°), whose 
particle spectrum does not significantly deviate from a 
thermal distribution. The suprathermal tail becomes 
flatter and more populated as the magnetic obliquity 
increases within the subluminal range, meaning that 



shocks close to the superluminality threshold 8 cr it « 34° 
accelerate more efficiently than quasi-parallel shocks. As 
8 increases from 0° to 30°, the slope of the power-law 
tail increases from —2.8 ± 0.1 to —2.3 ± 0.1 (subpanel 
a), the fraction of particles in the tail grows from ~ 1% 
to ~ 2% (subpanel b) and the energy fraction they ac- 
count for increases from ~ 4% to ~ 12% (subpanel c); 
at the end of our simulations the spectrum extends up 
to 7 m a X ~ 500 for 8 = 0° and reaches 7 max 2000 for 
8 = 30°. The negligible suprathermal tail observed for 
8 = 45° is much steeper than for subluminal shocks; as a 
conservative upper limit, the fraction of particles and en- 
ergy it accounts for is smaller by two orders of magnitude 
than in subluminal configurations. 

The dependence of the acceleration efficiency on mag- 
netic obliquity should affect the temperature A7 MB of the 
low-energy Maxwellian; we expect that, the larger the 
fraction of energy stored in the suprathermal tail, the 
cooler the low-energy bump, if the mean downstream 
particle Lorentz factor 7 were constant with magnetic 
obliquity. However, with increasing 8, a larger fraction 
of the upstream kinetic energy is converted to down- 
stream magnetic energy (see Table [TJ , so that the aver- 
age downstream particle energy monotonically decreases 
with 8. In summary, we expect that, within the sublumi- 
nal range, A7 MB should systematically decrease with 8, 
both because 7 decreases and because a larger fraction of 
downstream particle energy is stored in the suprathermal 
tail; across the boundary between subluminal and super- 
luminal configurations, A7 MB should sharply increase, if 
the drop in acceleration efficiency prevails over the de- 
crease of 7 with 8. Indeed, the thermal spread of the 
fitting Maxwellian decreases from ~ 4.6 for 8 = 0° to 
~ 4.2 for 8 = 30°, and then jumps back to ~ 4.5 for 
in Fig. 



8 = 45° 
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such trend can be seen in the posi- 
tion of the peak of the low-energy bump. 

In summary, the high-energy tail becomes harder 
and more populated as the magnetic obliquity increases 
within the subl uminal range, in agreement with semi 



analytica l (e.g., Kirk & Heavens 1989; Ballar d fc Heav- 
|ens||199l[) and Monte-Carlo (e.g . , [Bed narz & Ost rowskl 
P98nNIemiec fc Ostrowski|2004| JEUison fe Double|2004J ) 
calculations. For superluminal shocks, we find that seli- 
generated turbulence is not sufficiently strong to trap 
a significant fraction of particles in the vicinity of the 
shock for efficient acceleration. As we sh ow in Sj5] the 
upstream oblique waves discussed in §3. 2 1 regulate injec- 
tion into the acceleration process for 8 < 8 cr it . Since no 
waves are present in the upstream medium of superlumi- 
nal shocks, due to the absence of the stream of returning 
particles that trigger their growth, a negligible fraction 
of the incoming particles enters the acceleration process, 
and the acceleration efficiency is correspondingly small. 
The suprathermal tail barely seen for 8 = 45° in Fig. [lO] 
is probably populated by particles accelerated during a 
single shock crossing from upstream to downstream, and 
then advected downstream. We further explore the tran- 
sition between subluminal and superluminal shocks with 
regards to downstream particle spectra in Appendix B, 
where we show that the drop in acceleration efficiency 
across 6 CI n ~ 34° is rather abrupt. 

Finally, we remind that in Fig. [10] we are comparing 
downstream spectra computed in a slab whose distance 
behind the shock is the same for all magnetic obliqui- 
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Fig. 10. — Downstream particle energy spectra at time ui p t = 9000 for different magnetic obliquities: 9 = 0° (blue), 15° (green), 30° 
(red) and 45° (black). All spectra are computed in a slab 200 c/oj p -wide centered 400c/oj p downstream from the shock. The spectrum is 
normalized to the total number of particles within the considered slab, apart from a constant multiplication factor; it is supplemented with 
Poissonian error bars. Subpanels: a) Power-law slope of the high-energy tail as a function of the obliquity angle 9. The tail is defined to 
start at the intersection between the fitting Maxwellian and the fitting power law. b) Fraction of particles (in percent) in the high-energy 
tail, c) Fraction of energy (in percent) stored in the high-energy tail. 



ties. This may be misleading, since the acceleration pro- 
cess in subluminal shocks is not in steady state yet (see 
[4.2 1, and downstream particle spectra at a given time 
may vary among different longitudinal locations (e.g., 
compare panels d and e in Fig. uh. The particle en- 
ergy spectrum in a slab at fixed aistance downstream 
from the shock is affected by the efficiency both in ac- 
celerating particles and in transmitting downstream the 
freshly- accelerated particles. Since the latter decreases 
with increasing obliquity if cross-field diffusion is of mi- 
nor importance, shocks of larger obliquities may appear, 
as evaluated from a slab at fixed distance behind the 
shock, to accelerate less efficiently just because particle 
diffusion downstream from such shocks is slower. By 
comparing particle spe ctra in a slab at constant post- 



0=0° 



shock distance (Fig. 10 1, we find that subluminal shocks 
with higher obliquities accelerate more efficiently; this 
will hold true a fortiori when taking into account the 
obliquity dependence of the downstream diffusion prob- 
ability. 

4.2. Time evolution 



Figs. 11 12 and |T3| follow the time evolution of down- 



stream particle spectra for our representative sample of 
subluminal shocks (8 = 0°, 15° and 30° respectively), 
whereas Fig.[l4]refers to the superluminal shock 9 — 45°. 
In the upper subpanel of each figure, we plot the time 




Fig. 11. — Time evolution of the particle energy spectrum for 
9 = 0°: ujpt = 2250 (black), uj p t = 4500 (blue), u) p t = 6750 (green), 
uipt = 9000 (red). All spectra are computed in a slab 200 c/wp-wide 
centered 400 c/u p downstream from the shock. Upper subpanel: 
electron (blue) and positron (red) maximum Lorentz factor within 
the slab where the spectrum is computed, as a function of time. 
Lower subpanel: electron (blue) and positron (red) spectra at time 
u! p t = 9000 with Poissonian error bars, supplemented by the best- 
fitting Maxwellians (black dashed lines) and power laws (black solid 
lines) . 
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61 = 15° 



61=45° 





Fig. 12. — Time evolution of the particle energy spectrum for 
= 15°. See the caption of Fig. 1 1 1 1 for details. 



Fig. 14. — Time evolution of t he particle energy spectrum for 
8 = 45°. See the caption of Fig. |11| for details. Here, we do not 
attempt to fit the suprathermal tail. 



61=30° 




Fig. 13. — Time evolution of the particle energy spectrum for 
= 30°. See the caption of Fig. 1 1 1 1 for details. 



evolution of the electron (blue) and positron (red) max- 
imum Lorentz factor in the downstream slab where the 
spectrum is computed. The lower subpanel shows the 
electron (blue) and positron (red) spectra at the final 
time of our simulations (uj p t — 9000), together with the 
best-fitting low-energy Maxwellians (black dashed lines) 
and high-energy power laws (black solid lines). 

For both subluminal and superluminal magnetic con- 
figurations, the low-energy thermal bump does not sig- 
nificantly evolve with time. Instead, the high-energy tail 
grows with time as more and more particles are shock- 
accelerated, but its evolution is much more significant for 
subluminal shocks tha n for superluminal configurations. 
For 9 — 45° (Fig. [l4|, the particle spectrum for early 
times (black line) is remarkably thermal, and the minor 
suprathermal tail visible at later times stops evolving af- 



ter ujpt = 6650 (compare green and red lines in Fig. 14 1 



The upper subpanel in Fig. [14] confirms that the max- 
imum particle Lorentz factor saturates at 7 max ~ 120 
after ui p t w 5000. We are then reasonably confident that 



for 9 — 45° our simulation captures the full time evo- 
lution of the downstream particle spectrum, so that the 
claimed inefficiency of particle acceleration in superlumi- 
nal shocks should still hold true with a longer simulation 
timespan. 

For subluminal shocks, the high-energy tail stretches 
and flattens with time. Its ti me ev olution is moderat e for 
low obliquities (6 = 0° in Fig. [TT] and 6 = 15° in Fig.pl 
and more dramatic for magnetic inclinations c lose to the 
superluminality threshold (9 — 30° in Fig. 131. For 
9 = 0° the maximum particle Lorentz factor increases 
from 7 m ax ~ 300 at uj p t = 2250 to j max ~ 500 at 



UJ n t 



9000, 



whereas for 9 = 30° it grows much faster, 
from 7„ lax w 300 to 7 max ~ 2000. The fraction of parti- 
cles and energy stored in the tail also increases during our 
simulations. Therefore, for subluminal obliquities (and 
especially for angles 8 < C rit, which show the most sig- 
nificant time evolution) , the shock acceleration process is 
not in steady state yet, and the values we quote for the 
fractional contribution by the suprathermal tail to the 
downstream particle and energy census should be taken 
as lower limits. On the other hand, for late times (com- 
pare the green and red lines in Figs. 11 12 13 ), the shape 



of the spectrum does not change significantly, apart from 
the linear increase with time of the upper energy cutoff. 
Since the power-law spectral index is always p > 2, this 
means that the values we quote for the acceleration ef- 
ficiency at ujpt = 9000 will not significantly increase at 
later times. 

Concerning the time evolution of the maximum parti- 
cle Lorentz factor, the upper subpanels of Figs. [12] and 



13 (for 0=15° and 30° respectively) show a roughly lin- 
ear increase with time, whereas for 9 — 0° it saturates at 



7max ~ 500 after tv p t « 6000 (upper subpanel in Fig. 111. 
By running simulations of parallel shocks in which the 
first backward jump of the injector happens at different 
times (we tried 30000, 40000 and 50000 timesteps), we 
observe that, the sooner the injector starts jumping, the 
quicker the acceleration process saturates. Each jump 
unphysically removes the returning particles which are 
farthest upstream. This has two consequences: first, we 
artificially prevent these particles from returning to the 
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shock and possibly contributing to the downstream spec- 
trum; also, since the stream of returning particles trig- 
gers the growth of the upstream oblique waves, which 
play a major role in the acceleration process (see <f5|, 
we are unphysically affecting the whole process of shock- 
acceleration. In summary, the spectral time-evolution 
for 9 — 0° saturates due to our jumping injector, and 
the acceleration efficiency and maximum Lorentz factor 
we report for 9 = 0° should be taken as lower limits. 
Instead, for higher obliquities (0 — 15° and 30°), no un- 
physical saturation is observed in the spectrum, since 
the returning particles remain closer to the shock, so that 
fewer particles are removed by our jumping injector. The 
jumping injector does not seem to cause other unphysical 
effects on particle acceleration. 

Finally, we point out that the difference between the 
downstream spectra of pos itro ns and electrons for 9 = 
30° (lower subpanel of Fig. 13 at ui p t = 9000), which is 
remincscent of an analogous difference observed at earlier 
times (Fig.[7]e, at to p t = 2250), is an artifact of the small 
transverse size of our simulations. For magnetic obliq- 
uities near 9 — 30°, the stream of returning particles, 
which triggers the growth of the upstream oblique waves, 
does not propagate very far from the shock. As a conse- 
quence, the upstream pattern of electrostatic waves does 
not move far enough from the shock to let the shock elec- 
tric potential fluctuate in time (Fig.[6]A; for uj p t — 2250). 
Particles of one charge will then enter the shock with en- 
ergies systematically higher than particles of the opposite 
charge, and the spectral difference between electrons and 
positrons will persist over time. For smaller obliquities, 
the stream of returning particles recedes from the shock 
at higher speeds and the time-average of the shock elec- 
tric potential vanishes, so that no systematic difference 
is observed between the downstream spectra of electrons 
and positrons for 9 = 0° and 9 = 15° (lower subpanels in 
Figs. 11 and 12 1. For 9 = 30°, we have tested that the 



downstream spectral difference disappears when the sim- 
ulation box becomes larger (more than 4096 transverse 
cells) than the characteristic transverse coherence length 
of the upstream electrostatic oscillations. Therefore, no 
spectral asymmetry between electrons and positrons is 
expected in any realistic astrophysical scenario. More- 
over, the sum of electron and positron spectra does not 
vary among runs with different box width, proving that 
the shock electric potential does not play any appreciable 
role in the overall acceleration process. 6 

5. PARTICLE ACCELERATION MECHANISMS 

We explore the mechanisms responsible for particle ac- 
celeration in subluminal shocks by tracing the trajec- 
tories of representative high-energy particles extracted 
from the simulations. We find that particle accelera- 
tion is mostly mediated by Diffusive Shock Acceleration 
(DSA) for quasi-parallel shocks {9 < 10°), but Shock- 
Drift Acceleration (SDA) is the main acceleration mech- 
anism for larger, yet still subluminal, magnetic obliqui- 
ties. 

The trajectories shown in Figs. 



15 16 and tt7\ (for 



9 = 0°,15° and 30° respectively) are characteristic of 

6 As we comment in fjsj this is no longer true for electron-ion 
shocks, where ion acceleration via the so called Shock-Surfing Ac- 
celeration mechanism depends on the shock electric potential. 



positrons that end up with the highest energies. In each 
figure, we show the time evolution of the particle Lorentz 
factor (panel a) and 4-velocities (panel d) in the simula- 
tion frame. In panel a, red stands for negative y-velocity 
of the selected positron, yellow for positive y-velocity. 
Particle positions in the wall frame as a function of time 
are shown in panel c. To clarify the interaction between 
particles and electromagnetic fields, in panel b we plot 
the particle x-position relative to the shock (with the 
same color coding as in panel a), superimposed over two- 
dimensional strips of E y stacked in time. For every piece 
of the particle trajectory that lasts 22.5 Wp 1 , we plot a 
5 c/wp-wide E y strip which is centered along y at the par- 
ticle transverse location after 11.25 w" 1 since the start of 
that portion of the trajectory. Each strip extends along 
the whole x dimension. The electric field is still mea- 
sured in the wall frame but shifted along x so that the 
shock appears stationary. We choose to plot E y because 
it is the dominant electric component of the upstream 
oblique modes, and therefore it should be the primary 
agent of particle energization even for a parallel shock, 
where no uniform background E y<u is present. 

In order to estimate the relative contribution of DSA 
and SDA to the acceleration process in oblique shocks 
(9 0°), we make use of the fact that particles are 
shock-drift accelerated while they drift parallel (or anti- 
parallel, according to their charge) to the background 
motional electric field. Since in the wall frame of our 
simulations the downstream motional electric field £y,d 
is much smaller than the upstream E yu (see last col- 
umn m Table [l|, we may assume that in the downstream 
medium the particle energy stays approximately constant 
and we can relate the energy gain A7 SDA due to SDA to 
the upstream drift Ay u experienced by the particle: 



A7si 



(1) 



where q is the particle charge and E y<u = —(3oB ZtU < 
the upstream background electric field. Since the right 
hand side of eq. ([I]) can be computed from particle orbits, 
as well as the overall energy gain A7, the SDA fractional 
contribution A7 SDA /A7 will be a diagnostic of the rel- 
ative importance of DSA and SDA for the acceleration 
process. In the extreme limit of scatter-free SDA we ex- 
pect A7 SDA /A7 ~ 1. 

5.1. 9 = 0° 

The orbit of a representative high-energy positron from 
the simulation of a parallel shock is plotted in Fig. [T~5[ 
For u p t < 5800, the particle is approaching the shock 
from upstream with velocity (3{t) ss — /3 x (see x(t), the 
blue fine in panel c). Since y(t) stays approximately 
constant before the particle encounters the shock, the 
E y strips stacked in panel b are all centered at the same 
y-location for u p t < 5800. Therefore, panel b shows 
for uj p t < 5800 how the E y component of the upstream 
oblique modes (Fig. [2](? for u p t — 2250) evolves in time 
within a fixed longitudinal slice. It is thus evident that 
such waves are moving towards the shock and their phase 
velocity in the wall frame is smaller than the speed /3q « 1 
of the incoming particles. The positron energy var ies 
with time even before it encounters the shock (Fig. 15 a), 
in response to oscillations in 7/? x (i) (blue line in Fig.|15|rf) 
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Fig. 15. — Trajectory of a representative high-energy positron for 8 = 0°, in the wall frame of the simulation: a) Particle Lorentz factor 
as a function of time; colors show when the particle y- velocity is positive (yellow) or negative (red), b) Longitudinal position relative to 
the shock as a function of time, with the same color coding as in panel a. The particle orbit is superimposed over two-dimensional strips 
of E y stacked in time. For every piece of the particle trajectory that lasts 22.5 Up 1 , we plot a 5c/oj p -wide E y strip which is centered along 

y at the particle transverse location after 11.25 ujp 1 since the start of that portion of the trajectory. Each strip extends along the whole 
x dimension. ES stands for (E y /\E y \) \E y \ a . E y is measured in the simulation frame, but it is shifted in space so that the shock appears 
stationary, c) Particle orbit as a function of time (x(t) blue, y(t) red and z(t) green), relative to the particle position at its first encounter 
with the shock, d) Particle dimensionless 4- velocities 7/3 x (i) (blue), ■yf3 y (t) (red) and 7/3 z (t) (green) as functions of time. 



induced by the electrostatic component of the upstream 
waves. 

From uj p t ss 5800 to uj p t sa 6600, the selected positron 
scatters between the shock and the upstream medium, 
and its energy is ra pidly boosted up to a Lorentz fac- 
tor 7 w 250 (Fig. 15 a). As the particle energy grows, its 
gyro-radius also proportionally increases, as seen in panel 
c for y(t) and z(t). A moderate overall energy gain oc- 
curs upstream from the shock between oj p t ss 6600 and 
ujpt w 7500, and finally at ui p t w 7500 the positron is 
transmitted downstream, where its energy remains con- 
stant. Most of the accelerated particles end up down- 
stream, but the small fraction that escapes upstream 
populates the returning stream which triggers the for- 
mation of the upstream oblique waves. 

As seen in Fig. 15 b, the acceleration process involves 
stochastic bouncings between the shock and the up- 
stream medium; magnetic turbulence in the shock tran- 
sition layer prevents the particle from being transmit- 
ted downstream, whereas the upstream scattering is pro- 
vided by the oblique waves discussed in |3.2| whose 
amplitude is largest close to the shock. In particular, 



Fig. [T5]6 suggests that the positron gains energy when- 
ever its y- velocity f3 y {i) (red if negative, yellow if posi- 
tive) has the same sign as the electromagnetic field E y 
(black if negative, white if positive), and decelerates if 
they have opposite signs. The converse will be true 
for electrons (not shown). Overall boosting of the se- 
lected positron results from favorable (3 y {t)E y > scat- 
terings being more frequent than encounters in which 
Py(t)E y < 0. 

The stochastic nature of such wave-particle interac- 
tions, which are ultimately responsible for particle ac- 
celeration, is a defining feature of DSA. So, particle 
energization in parallel shocks proceeds via DSA, at 
the very least because the absence of a uniform back- 
ground motional electric field rules out SDA. However, 
in the traditional picture of first-order Fermi acceler- 
ation, accelerated particles scatter off magnetic turbu- 
lence embedded in the upstream and downstream media. 
Here, instead, high-energy particles bounce between the 
upstream medium and the shock front; moreover, the 
oblique modes responsible for the upstream scattering 
are not convected with the upstream plasma, as discussed 
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above. Nevertheless, the stochastic character of the ac- 
celeration process points towards a diffusive Fermi-like 
acceleration mechanism. Indeed, as theoretically pre- 
dicted for di ffusive acceleration in rel ativistic magnetized 
shocks (e.g., |Achterberg et al. 2001), we find that each 
acceleration cycle between the shock and the upstream 
medium involves an energy gain A7 ~ 7, at least for the 
limited number of shock crossings (< 10) observed in our 
simulation. 

These results for a parallel shock should also apply 
to quasi-parallel magnetic configurations, in which the 
obliquity is small enough that the upstream motional 
electric field does not si gnifi cantly affect the acceleration 
process. However, in §5.2| we show that SDA provides 
most of the energization for obliquity angles as low as 9 = 
15° , so that DSA is the dominant acceleration mechanism 
only for 9 < 10°. 

5.2. 9 = 15° 



Fig. 16 presents the acceleration process for a represen- 
tative positron in a shock with 9 = 15°. The particle is 
reflected back upstream on its first interaction with the 
shock (u! p t w 5100) and then gains energy at each shock 
encounter (ui p t 5200, uipt 5900 and uj p t « 8100). 
The selected positron eventually escapes upstream at 
time ujpt w 8200 with Lorentz factor 7 400, contribut- 
ing to the beam of returning high-energy particles. How- 



ever, due to the limited timespan of our simulation, we 
cannot exclude the possibility that the positron will fi- 
nally be deflected back and advected downstream. 

The time evolution of the positron Lorentz factor plot- 
ted in panel a suggests that energy gain is always as- 
sociated with f3 y (t) < (red), i.e., when the positron 
y-vclocity has the same sign as the uniform background 
electric field E y , u < 0. Whenever the positron gyro- 
orbit is fully contained in the upstream region (e.g., from 

6000 to 



ujpt ps 5300 to ujpt 



5800, or from ujpt 



ujpt w 8000, as well as after the selected positron es- 
capes upstream), the energy gain within the "red" half 
of the gyro-period is balanced by an approximately equal 
energy loss during the "yellow" half, since the back- 
ground motional electric field -By jU stays approximately 
constant throughout the gyro-orbit. However, when the 
positron trajectory intersects the shock front, in the "up- 
stream" portion of its clockwise Larmor gyration around 
the oblique field the positron y-velocity has always the 
same sign as the background motional electric field, thus 
resulting in systematic acceleration. Instead, the "down- 
stream" portion of the orbit does not result in decelera- 
tion because in the simulation frame the motional elec- 
tric field downstream is much smaller than in the up- 
stream. 7 For electrons (not shown), whose Larmor gyra- 
tion around the magnetic field is counter-clockwise, the 

7 We emphasize that a similar line of reasoning can be made 
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y-velocity in the upstream medium will have opposite 
sign than -E y ,u, resulting again in net energization. 

In summary, particles are systematically accelerated 
by the upstream motional electric field whenever their 
gyro-orbit intersects the shock surface, as expected for 
the Shock-Drift Acceleration process. SDA is inevitably 
associated with a VB drift along (or opposite, according 
to the particle charge) the background electric field, due 
to the magnetic field jump at the shock. This is seen, 
for the positron in Fig. |16| in the slow drift along — y of 
y(t) in panel c. Such a drift will be more significant for 
higher obliquities, since the magnetic field jump at the 
shock is larger. 

A comparison between x(t) and z(t) in Fig. 16 c sug- 
gests that our representative positron is heading up- 
stream along the oblique magnetic field. As compared 
to a parallel shock, it will take longer to escape ahead 
of the shock, the more so for higher obliquities. Many 
gyro-orbits could then intersect the shock during a single 
shock-crossing of the positron gyro-center, thus result- 
ing in multiple energizations (as happens at uj p t rs 5200 
and ojpt sa 5900). The particles that escape upstream 
are not systematically shock-drift accelerated any longer, 
until they return to the shock (w p t ~ 8000 in Fig. 16 1, 



in any other frame. In the particular frame of the simulation, the 
downstream motional electric field is E y ^ <C Ey,u, which simplifies 
the argument. However, the shock is moving in this frame. 



due to deflection by the upstream background fields or 
scattering off the upstream oblique waves. In that re- 
spect, even though the upstream motional field E y u is 
the primary agent for particle energization, DSA coop- 
erates with SDA to trap the accelerated particles at the 
shock, thus increasing the overall acceleration efficiency. 
Making use of eq. (IT| with a sample of ~ 1000 high- 
energy particles, we are able to confirm that the ratio 
A7 SDA /A7 is in dependent of the total energy gain A7 
(as predicted by |Jokipii 1982), and we estimate the frac- 
tional contribution by SDA to be ~ 60%. Therefore, for 
obliquity angles as small as 9 — 15°, SDA is already the 
main acceleration mechanism. 



5.3. 6» = 30° 

The acceleration process for a representative positron 
in a shock with 9 — 30° is shown in Fig. [IT] As the obliq- 
uity angle approaches the critical boundary between sub- 
luminal and super luminal configurations (# cr it ~ 34°), 
particles reflected upstream can hardly escape ahead of 
the shock. They are effectively trapped at the shock 
while heading upstream along the oblique magnetic field, 
since for 9 < 6* cr it their gyro-center a;- velocity is compara- 
ble to the shock speed. Before escaping upstream, their 
orbit may intersect the shock many times (many more 
than for 9 — 15°), thus resulting in repeated energiza- 
tions via SDA. 
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As discussed in £5.2 for a particle gyrating across the 
shock any excursion into the upstream region results in 
efficient SDA by the upstream motional electric field E YyU 
(compare panels a and b in Fig. 17 1, whereas in the down- 



stream medium its energy stays roughly constant since 
E y $ <C E y u in the simulation frame. The correspond- 
ing time evolution of the Lorentz factor 7(2) resembles 
a stairway (see Fig. 17 a), with steps that grow linearly 
with 7 since the time spent upstream is proportional to 
the particle energy. When a particle completes a full 
gyro-orbit in the upstream medium (betw een tu p t ~ 7200 
and ujpt w 7600 for the positron in Fig. [l7| , its energy 
gain is balanced by an equal energy loss. Due to scat- 
tering off the incoming oblique waves or regular deflec- 
tion by the upstream background fields, the particle may 
be re-captured by the shock, where its acceleration via 
SDA resumes. The representative positron in Fig. [17] 
as most of the high-energy particles extracted from our 
9 = 30° simulation, is still being shock-accelerated at the 
final simulation time uj p t — 9000; its final Lorentz factor 
7 w 1000, which is already much larger than in shocks 
with smaller obliquities, is therefore just a lower limit. 

Using eq. (fTl) with a sample of ~ 1000 high-energy par- 
ticles, we fina that on average SDA contributes a fraction 
~ 90% to the overall energy gain, so that diffusive scat- 
tering plays a minor role in the acceleration process for 
9 = 30°. Thus, with increasing magnetic obliquity in 
the subluminal range, the acceleration process switches 
from DSA, which is the only mechanism acting for par- 
allel shocks, to SDA, which for 9 < 9 crlt operates in an 
essentially scatter-free regime. 

For 9 = 30° the upstream oblique waves are important 
to regulate injection into the acceleration process: all the 
particles that end up gaining the highest energies have 
been deflected from their upstr eam straight path before 
encountering the shock (see Fig.|17|fe at u> p t 6500). We 
have tested that incoming particles enter the acceleration 
process only if their first shock encounter results in a 
reflection back upstream, whereas transmitted particles 
populate the low-e nergy bump of downstream spectra. 
In agreement with |Ballard fc Heavens ( |1991[ ) , we have 
verified that such reflected particles had encountered the 
shock with p itch angle cosine /x H t in the de Hoffmann- 
Teller frame ( de Hoffmann fc Teller|1950 l which satisfies 
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where i?d,HT / B UjHT is the magnetic field jump in the de 
Hoffmann- Teller frame. In other words, only the incom- 
ing particles that are deflected by the upstream waves, 
such that their pitch angle satisfies eq. p I when encoun- 
tering the shock, will be reflected back upstream and 
participate in the SDA process. 

In summary, while for quasi-parallel shocks (9 < 10°) 
the upstream oblique waves generated by the beam of 
returning particles provide the upstream scattering for 
DSA, for magnetic obliquities close to 6* cr ; t the same 
waves mediate the injection process into SDA. 

6. DISCUSSION AND SUMMARY 

We have explored by means of 2.5D PIC simulations 
the internal structure and acceleration properties of rel- 
ativistic collisionless shocks propagating in a magnetized 



pair plasma. A strong magnetic field may significantly 
suppress cross-field motion an d therefore quench the 
growth of the Weibel instabil ity (Wcibel 1959 Mcdvedev 
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pendicular shocks are triggered by magnetic reflection of 



the in coming particles from th e shock-compressed fields 
(e.g., lAlsop fc Arons|[l988l ). |Spitkovsky| p005l finds 
that the transition between weakly magnetized (Weibel- 
mediated) and strongly magnetized (reflection-mediated) 
perpendicular shocks happens around an upstream mag- 
netization a ~ 10~ 3 . In this work, we focus on strongly 
magnetized upstream flows [a — 0.1) and we study the 
dependence of the shock properties on the inclination an- 
gle 9 between the upstream magnetic field and the shock 
normal, as measured in the "wall" frame of our simula- 
tions. 

In magnetized shocks, particle gyro-centers are con- 
strained to move along the magnetic field. It is then 
natural to define a critical obliquity angle 9 cl a such that 
in "superluminal" shocks (9 > C rit) particles cannot es- 
cape upstream along the magnetic field lines. If the up- 
stream plasma has Lorentz factor 70 = 15 and magne- 
tization a = 0.1, MHD calculations yield # cr it ~ 34° 
in the simulation frame, weakly dependent on both 70 
and a. Since for 9 > 6 cr - lt particles are quickly ad- 
vected downstream from the shock, the acceleration ef- 
ficiency should be strongly suppressed for superluminal 
configurations. In fact, we find that magnetized rela- 
tivistic pair shocks produce a significant population of 
suprathermal particles only for "subluminal" magnetic 
geometries (9 < 9 cl . it ). This holds true for a wide range 
of upstream Lorentz factors (we also tried 70 = 5 and 
70 = 50, for fixed a — 0.1) and magnetizations (we also 
tested a = 0.03 and a = 0.3, for fixed 70 = 15). We con- 
clude that only a narrow range of obliquity angles allows 
for efficient particle acceleration in relativistic (70 > 5) 
magnetized [a > 0.03) pair shocks: as measured from 
the upstream frame, the magnetic field should lie within 
a cone of half-opening angle 34° /70 around the shock 
normal to produce a significant population of nonthermal 
particles. 8 

For quasi-parallel magnetic configurations (9 < 10°), 
some of the incoming particles are scattered back up- 
stream and then gain energy by bouncing between the 
shock and the upstream medium. The upstream scat- 
tering is provided by oblique waves which propagate to- 
wards the shock and are self-consistently triggered by the 
shock-accelerated particles that escape upstream; mag- 
netic turbulence in the shock transition layer reflects the 
high-energy particles back upstream. The stochastic na- 
ture of the acceleration process suggests that DSA is the 
main energization mechanism in quasi-parallel magne- 
tized pair shocks, at the very least because Shock-Surfing 
Acceleration does not operate in electron-positron shocks 

8 We remark that this conclusion applies only to electron- 
positron shocks. For electron-ion magnetized shocks, the electric 
potential across the shock complicates the particle dynamics and 
we cannot exclude that shocks which are nominally superluminal 
may be moderate particle accelerators. Relativistic shocks prop- 
agating in magnetized electron-ion plasmas will be discussed in a 
forthcoming paper. 
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and the absence of a background motional electric field 
rules out SDA in parallel shocks. For 9 — 0°, we find 
that at the final time of our simulation (uj p t — 9000) the 
suprathermal tail of the downstream energy spectrum 
accounts for ~ 1% of the downstream particle number 
and ~ 4% of the energy, with a power-law spectral in- 
dex — 2.8±0.1; however, these may be lower limits, since 
further evolution of the spectrum in our simulations of 
quasi-parallel shocks may be artificially suppressed by 
our choice of boundary conditions. 

Whereas for quasi-parallel magnetic geometries (9 < 
10°) the shock-generated turbulent fields are typically 
stronger than the uniform motional electric field, and so 
SDA plays a minor role in the acceleration process, its 
mean fractional contribution to the particle overall en- 
ergy gain increases to ~ 60% for 9 = 15° and reaches 
~ 90% for 9 — 30°. Therefore, the acceleration process 
for 9 < 6 CT it is almost scatter-free: particles heading up- 
stream along the oblique magnetic field hardly escape 
from the shock, and they can experience multiple shock- 
drift energizations as their orbit repeatedly crosses the 
shock surface. However, we cannot exclude that a frac- 
tion of the shock-drift accelerated particles may even- 
tually gain enough energy to become effectively "un- 
magnetizcd" , i.e., their gyrocenter would no longer be 
constrained to follow the background field. Such parti- 
cles may be further accelerated by a diffusive Fermi-like 
mechanism, and the SDA process would then serve as a 
pre-acceleration stage before the onset of DSA. The up- 
stream oblique waves, which provide the upstream scat- 
tering in quasi-parallel shocks, may also mediate this fur- 
ther step of diffusive acceleration for larger, but still sub- 
luminal, obliquities. We also find that in 9 < 9 crlt shocks 
these waves regulate the process of particle injection into 
SDA. 



As predicted by Jokipii 

< 



([1982], we find that SDA in 
y cr it shocks is a much faster process than DSA in 
quasi-parallel shocks, since each SDA cycle takes approx- 
imately one gyro-period, which is much shorter than any 
diffusive scattering timescale. As the acceleration mech- 
anism switches from DSA to SDA with increasing obliq- 
uity, the acceleration efficiency increases as well. For 
9 = 30°, the number and energy fractions contributed 
by the suprathermal tail of the downstream particle spec- 
trum are respectively ~ 2% and ~ 12% and its power- 
law spectral index is —2.3 ± 0.1. Since acceleration in 
9 k, # cr it shocks is close to a steady state around the fi- 
nal time of our simulations (ujpt = 9000), the efficiencies 
quoted above should be near the saturation values. 

For superluminal shocks (9 > C rit)j we find that if 
the upstream background field is not highly turbulent on 
scales smaller than the transverse size of our simulations 
(~ 100c/cjp), the self-generated turbulence is not strong 
enough to trap a significant fraction of particles at the 
shock for efficient acceleration; rather, most of the parti- 
cles are quickly advected downstream with the magnetic 
field lines. Monte-Carlo simulations of relativistic mag- 
netized shocks found eviden ce for particle acceleration in 



perpendicular shocks (e.g., |Bednarz fc Ostrowski 



Niemiec fc Ostrowskif2 004 ; Elli son fc Doubl e 200 



ki 
34 



1998 
, but 



they required ex ceptionally high levels of turbulence ( Os 
|trowski fc Bednarz||2002[ ) . We show via PIC simulations 
with laminar upstream fields that superluminal shocks 
do not self-consistently generate such large degrees of 



turbulence, and they are inefficient particle accelerators. 
In Appendix B we show that the drop in acceleration 
efficiency at # cr it ~ 34° between subluminal and super- 
luminal shocks is rather abrupt. 

We have tested that our results do not significantly 
depend on the number of particles per computational 
cell (we tried 2 and 8 particles per species per cell) or 
the transverse size of our two-dimensional simulation box 
(we tried 512, 1024, 2048 and 4096 transverse cells). We 
have verified that the geometry of the magnetic field with 
respect to the simulation plane (in-plane or out-of-plane 
magnetic field) does not change our main conclusions (see 
Appendix A); moreover, our 2.5D results have been con- 
firmed by running a limited suite of simulations in three- 
dimensional domains with 64, 128 or 256 cells along each 
transverse dimension. 

A direct comparison of 2.5D simulations with 3D runs 
is possible, since for moderate obliquities (9 < 45°) the 
gyro-orbits of downstream particles mostly lie outside 
the two-dimensional simulation plane, and the flow ther- 
malizes behind the shock to a three-dimensional distri- 
bution. For the obliquity angles investigated with 3D 
simulations (9 = 30° and 9 — 35°), we see a remark- 
able agreement with the results of our 2.5D runs (for 
9 = 30°, see Appendix A). |Jones et al.| ( |1998| found that 
numerical simulations in two spatial dimensions could ar- 
tificially inhibit particle diffusion across field lines, which 
is required for efficient particle acceleration in superlu- 
minal shocks. Instead, we confirm with 3D simulations 
that superluminal relativistic shocks propagating in pair 
plasmas with laminar fields do not produce a significant 
population of nonthermal particles, at least for the do- 
main sizes we tried in 3D. 

In summary, our results underscore that particle ac- 
celeration in relativistic (70 > 5) magnetized (a > 0.03) 
pair flows is efficient only for magnetic obliquities 9 < 
9 crit w 34°. As seen from the upstream fluid frame, only 
nearly-parallel shocks (with obliquity < 34° /70) produce 
a significant suprathermal tail of accelerated particles. 
These findings place constraints on the models of Pulsar 
Wind Nebulae (PWNe), Active Galactic Nuclei (AGN) 
jets and Gamma-Ray Bursts (GRBs) that invoke particle 
acceleration in relativistic magnetized shocks. 



sarov & Lyubarsky 


2003 


2004 


have shown that polar jets wit 



MOT 



Del Zanna et aJ 
1 the observed velocities 
can develop provided that the magnetization of the rela- 
tivistic (70 > 10 6 ) electron-positron pulsar wind is large 
enough (a > 0.01) to collimate the jet by hoop stresses 
downstream of the wind termination shock. Since the 
asymptotic wind magnetic field is mostly toroidal around 
the pulsar axis, the termination shock will be quasi- 
perpendicular. If the wind magnetization is a > 0.03, 
our results predict that the termination shock should not 
produce a significant population of nonthermal particles, 
which is required instead by the observed synchrotron 
and inverse Compton emission. However, for oblique ro- 
tators, the pulsar wind in the equatorial plane should 
consist of stripes of alternating field, and the field may 
dissipate via magnetic reconnection of the opposite polar- 
ity stripes before the ter mination shock (|Coroniti||1990| 



Lyubarsky & Kirk 20011; the effective magnetization of 
the shock would then be lower and particle acceleration 
would proceed like in an essentially unmagnctizcd shock. 
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Alternatively, magnetic reconnection of the striped wind 
field when compressed at the termin ation shock may be 



responsible for particle energizat ion (Petri & Lyubarsky 



2007| |Lyubarsky fc Liverts||2008[ 

Synchrotron and inverse Compton emission from the 
core of blazar jets is usually attributed to high-energy 
electrons and positrons accelerated in mildly-relativistic 
internal shocks (70 ~ 2 in the jet comoving frame). 
Such shocks should be quasi-perpendic ular, since po- 



larization measurements of radio knots (Gabuzda et al. 
|2004| |Pushkarev et al.||2005[ ) indicate that the magnetic 
field is mostly transverse to the jet axis. If the jet flow 
is significantly magnetized (a > 0.03), the presence of 
a substantial population of high-energy particles acceler- 
ated at such perpendicular magnetized shocks would be 
hard to explain, unless particle acceleration proceeds not 
in the central "spine" of the jet but rather at its edges; 
here, interaction of the jet with the surrounding medium 
may create subluminal configurations where particle ac- 
celeration would be more efficient. Otherwise, the jet 
flow should be weakly magnetized (a < 0.03). Alterna- 
tively, particle acceleration may occur not in shocks but 
in reco nnection layer s of a Poynting-dominated jet, al- 
though Sikora et al. ( 2005 1 have shown that the jet is 
likely to be matter-dominated at the pc-scale emission 
region. 

Similarly, in the standard "fireball" scenario of GRBs, 
most models of the prompt emission require a population 
of energetic particles accelerated at mildly-relativistic in- 
ternal shocks. If the ejecta are strongly magnetized, pos- 



sibly by the field of the progenitor star, and the magnetic 
field is mostly transverse to the outflow, we would predict 
inefficient particle acceleration and thus negligible non- 
thermal emission, contrary to observations. This could 
be circumvented if the relativistic outflow primarily takes 
the form of a large-scale Poynting flux, and particle ac- 
celeration occurs not in shocks but by dissipa tion of mag- 
netic energy via current-driven instabilities (Lyutikov & 
Blandford 2003|. 



In summary, our results suggest that models of non- 
thermal emission in PWNe, AGN jets and GRBs which 
require particle acceleration in relativistic pair shocks 
should resort to either nearly-parallel or weakly mag- 
netized shocks. Alternatively, other acceleration mecha- 
nisms, e.g., by dissipation of the field via magnetic recon- 
nection, should be invoked, or the upstr eam flow should 
conta i n a significant ionic c omponent (Hoshino et al. 
1992| Amato & Arons 2006) or small-scale turbulence 
( Sironi fc Goodman||2007| ). 
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APPENDIX 

A) 2.5D SIMULATIONS WITH IN-PLANE UPSTREAM BACKGROUND MAGNETIC FIELD 

In the 2.5D simulations presented in this work, the fact that the upstream background magnetic field lies in a plane 
which is perpendicular to the simulation plane ( "out-of-plane" magnetic configurations) may be suspected in reducing 
turbulence along the field. Here we test this by analyzing the results of 2.5D simulations with an "in-plane" upstream 
blackground magnetic field, for each of the subluminal oblique configurations presented in Sj3j i.e., 9 — 15° and 9 = 30°. 

In the in-plane case, the upstream background magnetic field will have components -B x ,u > and B y u > (so that 
9 = arctan[fly jU /S XjU ]), and the associated upstream motional electric field will be E u — /3o-B y , u z. As regards to the 
electromagnetic fields, the direction +z in out-of-plane configurations should then correspond to in- plane +y, and 
out-of-plane +y to in-plane — z. Apart from this overall rotation of the coordinate system, we expect out-of-plane 
and in-plane geometries with the same magnetic obliquity to show comparable results in terms of shock structure and 
particle spectra. Any difference should be an artifact of the reduced dimensionality of our simulation box and would 
be a caveat against the generalization of our 2.5D results to three-dimensional problems. 



Shock structure 



Figs. |18| and fl9| show the shock structure and positron phase space at time uj p t — 2250 for in-plane magnetic^ 



configurations with obliquity 9 — 15° and 9 = 30° respectively. Their out-of-plane counterparts are plotted in Figs. [3] 
and [6] respectively. We find that the shock internal structure does not significantly change by varying the orientation 
of the field with respect to the simulation plane. 
Out-of-plane and in-plane configurations with the same obliquity show a comparable shock velocity and similar 



= 15° CJ„t=2250 
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1000 



Fig. 18. — Shock structure and positron phase space at time ui p t = 2250 for 8 = 15°, with upstream background magnetic field lying in 
the simulation plane. The fluid quantities are normalized as in Fig. [2] but here B z and E y are in units of the transverse component -By jU 
of the upstream background magnetic field. 
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Fig. 19. — Shock structure and positron phase space at time u p t = 2 250 for = 30°, with upstream background magnetic field lying in 
the simulation plane. The fluid quantities are normalized as in Fig. |18| 



transversely-averaged density (panel b) and magnetic energy (panel d) profiles, as relates to the width of the shock 
transition region (~ 10c/o-> p ), the jump across the shock (in agreeme nt with three-dimensional MHD calculations) 



and Fig. [19]for 9 = 30°). For 9 = 15°, the y-averaged 
le upstream kinetic energy density, both for out-of-plane 



and the presence of a fast and a slow shock (compare Fig. k3 
magnetic energy profile peaks in the shock layer at ~ 30% of t 

(Fig. |4]d) and in-plane (Fig. 18 d) configurations, and in both cases it decays downstream on a lenghtscale ~ 50c/w p . 
For 9 = 30°, the overall fluid structure does not significantly change between in-plane and out-of-plane configurations, 
but a few minor differences ar e p resent: for the in-plane case, the magnetic energy in the shock layer slightly exceeds 
the downstream plateau (Fig. Il9]rf), which was not observed in its out-of-plane counterpart (Fig. p\d); also, the speed 
of the slow shock, tracked from the point where the magnetic energy starts to decrease downstream from the main 
shock, seems to be higher for in-plane than out-of-plane co nfig urations. For both 9 = 15° and 9 
electric field (panel h in Figs. B] and [6| not shown in Figs, 
due to a residual downstream Dulk velocity transverse to t 



30°, the motional 
18 and 19 1 does not drop to zero behind the main shock 
re upstream flow {P Yl d < and (3 Z 



; d < for in-plane and 

out-of-plane configurations respectively) . The transient electromagnetic precursor wave generated by coherent Larmor 
bunching at the shock, which showed up in B z and E y for out-of-plane oblique configurations, is now seen in the 



corr esp onding in-plane fields B y and E z with comparable wavelengths and amplitude (further upstream than Figs. 18 
and 19 show). 

The phase space plots in panels l-n of Figs. 18 and 19 are similar to their counterparts in Figs. Hand [6] respectively, 
provided that out-of-plane 7/3 z and 7/3 y are compared with in-plane 7/3 y and — 7/3 z respectively. As observed for out- 
of-plane configurations, the returning high-energy particles slide along the upstream ma gnetic fi eld, which means for 
the in-plane case that they show large and positive 7/3 x and 7/3 y (panels I and m in Fi gs. 18 19|. The in-plane excess 

is due to V-B drift, 



of positrons with 7/? z > in the upstream region close to the shock (panel n in Figs. 18 



ancTT 
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(a) For fixed magnetic obliquity 9 = 15°, downstream particle 
spectra at time u) p t = 9000 for different magnetic field geometries: 
magnetic field lying either in the simulation plane (black) or in a 
plane perpendicular to the simulation plane (red). 



(b) For fixed magnetic obliquity 8 = 30°, downstream particle 
spectra at time io p t = 9000 for different magnetic field geometries: 
two-dimensional simulations with either in-plane (black) or out-of- 
plane (red) magnetic field and a three-dimensional run (green). 



Fig. 20. — Comparison of downstream particle spectra among different geometries of the upstream background magnetic field, for 
magnetic obliquity 9 = 15° in panel (a) and 9 = 30° in panel (b). 

and it corresponds to the asymmetry towards negative values of 7/3 y (panel m in Figs. [4] and [6] respectively) observed 
for out-of-plane configurations. 

Aside from the similarities between in-plane and out-of-plane results, Figs. [18] and [l9| improve our knowledge of the 
internal structure of oblique shocks. The two-dimensional plots of magnetic energy density (panel c) and B z (panel 
e) of Figs. [4] and [6] show islands of low magnetic field in the downstream region of the strong shock. Fewer islands are 
produced at the shock at later times; when formed, each island slowly decays with time at constant x-location, and it 
disappears when caught up by the slow shock 



18] and fl9| clarify that such islands 
as confirmed also by 3D 



The corresponding panels in Figs 
are actually tubes of low magnetic energy stretching along the downstream magnetic fieL 
simulations with 256 cells along each transverse dimension. 

Most importantly, the electromagnetic component of the upstream oblique waves, which was mostly seen in B z and 
E y for the out-of-plane case (panels e, g in Figs.j4j6|, shows up mainly in the same field components also for in-plane 
configurations. In other words, the wave magnetic field is still mostly perpendicular to the simulation plane, regardless 
of the orientation of the upstream background magnetic field with respect to the simulation plane. A further study of 
such modes will therefore have to be fully three-dimensional to capture their true physical nature. What Figs. [18| and 
19 (panels e, g, i) add to our understanding of such waves is that, in the region far upstream where they are triggered, 
the surface of constant phase is roughly a cone of constant opening angle around the upstream background field. 

Particle energy spectra 



Fig. 20 shows at Lu p t — 9000 a comparison of downstream particle spectra between in-plane (black line) and out- 



of-planc(red line) configurations: for 9 = 15° (panel a), the spectrum is remarkably similar between the two cases, 
whereas a significant difference is seen for 9 — 30° (panel 6), with the high-energy tail for the out-of-plane configuration 
being much more pronounced than in its in-plane counterpart (and the peak of the low-energy Maxwellian being 
correspondingly cooler). However, the out-of-plane spectrum for 9 = 30° is in g ood agreement with a 3D simulation 
having 64 cells along each transverse dimension (green line in panel b of Fig. 20 1, which suggests that results from our 



2.5D out-of-plane simulations may be confidently applied to realistic three-dimensional scenarios. 

Concerning the fact that the 9 = 30° spectrum at uj p t — 9000 shows a significant difference at high energies between 
in-plane and out-of-plane configurations, we report that the spectrum computed in the shock transition layer (not 
shown) is remarkably similar between the two cases, suggesting that it is not the intrinsic acceleration efficiency to 
differ but rather the transmission probability of shock-accelerated particles into the downstream medium. 

B) THE SHARP TRANSITION IN ACCELERATION EFFICIENCY BETWEEN SUBLUMINAL AND SUPERLUMINAL 

SHOCKS 

In Fig. [TO] we have shown that only subluminal shocks present a pronounced suprathermal tail in downstream particle 
spectra, whereas the particle distribution function behind superluminal shocks does not significantly deviate from a 
Maxwellian. Here, we show that the drop in acceleration efficiency between subluminal and superluminal configurations 
happens remarkably close to the boundary 9 cr it « 34° predicted by MHD calculations, and it is rather abrupt. 

In Fig. [21] downstream spectra for 9 = 28° (violet), 30° (red), 31° (black), 32° (blue), 35° (green) and 45° (yellow) 
are compared for different times. Within the subluminal range, lower obliquities (9 = 28°) show the most significant 
suprathermal tail at early times, whereas shocks with magnetic inclinations closer to 9 cr a (6 — 30° and 31°) are the 
most efficient accelerators at later times. The closer to the superluminality threshold, the faster the time evolution of 
the spectrum, in agreement with the fact t hat particle acceleration in shocks with 1 



< 



'crit 



is mostly mediated by SDA, 

which is a very fast process ( Jokipii 1982). More importantly, the onset of efficient SDA is extremely rapid, since it is 
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driven by a positive feedback, as discussed in [5.3 as the number of accelerated particles which escape upstream gets 
large enough to trigger the generation of the upstream oblique waves, a larger fraction of the incoming particles will 
be deflected from their straight path and their pitch angle at the shock will more likely satisfy the relation in eq. ([2 1. 
More particles will then be reflected upstream and injected into the acceleration process, in a runaway cascade towards 
larger acceleration efficiency The dramatic difference observed in the high-energy tail of 9 = 31° between oj p t = 2250 



and 



wj, = 4500 and of 9 = 



32° between ui p t = 6750 and w p t — 13500 (thin blue line in the bottom right panel of 



Fig. 21 ) corresponds to the onset of efficient SDA. In response to the increased acceleration efficiency, the peak of the 



low-energy Maxwellian shifts to lower temperatures. 

Since the positive feedback required for the onset of efficient SDA relies on upstream waves triggered by the returning 
particles, we expect that in superluminal shocks, where particles are not able to return upstream along the magnetic 
field, acceleration should be strongly suppressed. Indeed, the downstream spectra of the superluminal shocks 9 — 35° 
and 6 = 45° are mostly thermal, although a minor high-energy tail, further suppressed as 9 increases, is observed at 
late times. Since spectra computed ahead of the shock (not plotted) do not show any sign of high-energy particles, 
we believe that the residual suprathermal tail observed downstream is probably populated by particles shock-drift 
accelerated during a single shock crossing from upstream to downstream and then advected downstream. We are 
confident that these conclusions hold true for times longer than the timespan of our simulations, for the following 
reasons: i) between m D t = 6750 and uj p t — 9000, the downstream spectra of superluminal shocks do not significantly 
evolve (see also Fig. 14), so that they are presumably close to a steady-state; ii) in the whole timespan explored by our 



simulations, there is no evidence for energetic particles escaping ahead of superluminal shocks, so that no upstream 
waves should be generated and the onset of efficient SDA should never occur. 

The subpanels in the last panel of Fig. 21 show at time uo p t = 9000 the power-law slope of the suprathermal tail 
and the fraction of particles and energy stored in the tail. We see that, with increasing magnetic obliquity within the 
subluminal range, the high-energy tail becomes flatter, more populated and more energetic. We remark that 9 = 32°, 
which is still subluminal, does not accelerate very efficiently yet at uj p t = 9000, but it will fi nall y show a suprathermal 
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for ujrd = 



13500). We 
= 9000 and 



tail similar to 9 — 30° and 9 = 31° (see the thin blue line in the bottom right panel of Fig. 
also report that the downstream spectra for 9 — 30° and 9 = 31° do not significantly evolve"T5etween uo p l 
uo p t — 13500 (not shown), apart from the linear increase of the high-energy cutoff. This suggests that the timespan of 
our simulations can reasonably capture the evolution of 9 < 9 CI i t shocks until they reach a steady-state, and that the 



values we quote in the subpanels of Fig. K21J (and Fig. 10 as well) are close to the saturation values 
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a +=4500 




Fig. 21. — Time evolution of downstream particle spectra for different obliquities, across the boundary cr ; t « 34° between subluminal 
and superluminal shocks: 9 = 28° (violet), 30° (red), 31° (black), 32° (blue), 35° (green) and 45° (yellow). The subpanels a)-c) in the last 
panel show at time w p t = 9000 the power-law slope of the suprathermal tail and the fraction of particles and energy stored in the tail, as 
a function of the obliquity angle 9. The thin blue line in the bottom right panel shows the particle spectrum for 9 = 32° at u) p t = 13500. 



